/*

bvqx_srf_parseneighborsc.c  - parse neighbors more quickly

FORMAT:       [neighbors, skip] = bvqx_srf_parseneighborsc(neid, numvtx)

Input fields:

      neid        neighbor data
      numvtx      number of vertices

Output fields:

      neighbors   Nx2 cell array with neighbors
      skip        number of elements to skip in input

% Version:  v0.7g
% Build:    9101521
% Date:     Oct-15 2009, 9:27 PM CEST
% Author:   Jochen Weber, SCAN Unit, Columbia University, NYC, NY, USA
% URL/Info: http://wiki.brainvoyager.com/BVQXtools

 */

#include "mex.h"
#include "matrix.h"

void mexFunction(int nlhs, mxArray *plhs[], int nrhs, const mxArray *prhs[])
{

    int nv = 0, ne = 0, skip = 0, nn, vn = 0, vc;
    const int *ind;
    double dnv = 0.0;
    mxArray *outa;
    double *outd;
    const mxArray *cell;
    const double *dind;
    
    /* argument check */
    if (nrhs != 2 || nlhs != 2)
		mexErrMsgTxt("Bad number of input/output arguments.");
    if (!mxIsUint32(*prhs) ||
        (mxGetClassID(prhs[1]) != mxDOUBLE_CLASS))
		mexErrMsgTxt("Input arguments must be of type uint32/double.");
    if (mxGetNumberOfElements(prhs[1]) != 1)
		mexErrMsgTxt("Number of vertices must be scalar.");
    dind = (const double*) mxGetPr(prhs[1]);
    if (mxIsInf(*dind) ||
        mxIsNaN(*dind) ||
        (*dind < 3.0))
		mexErrMsgTxt("Invalid number of vertices.");
    skip = nv = (int) *dind;
    ne = mxGetNumberOfElements(*prhs);
    if (ne < (4 * nv))
		mexErrMsgTxt("Neighbors stream too short.");
    ind = (const int*) mxGetPr(*prhs);
    
    /* prepare output */
    *plhs = mxCreateCellMatrix(nv, 2);
    if (*plhs == NULL)
        mexErrMsgTxt("Error creating Neighbors cell array.");
    plhs[1] = mxCreateDoubleMatrix(1, 1, mxREAL);
    if (plhs[1] == NULL)
        mexErrMsgTxt("Error creating skip indicator.");
    if (((double*) mxGetPr(plhs[1])) == NULL)
        mexErrMsgTxt("Error allocating skip indicator memory.");
    
    /* loop over vertices */
    for (vc = nv; vc > 0; --vc) {
        
        /* get number of neighbors for vertex */
        nn = *ind++;
        skip += nn;
        
        /* check skip vs ne ! */
        if (skip >= ne)
            mexErrMsgTxt("Too little data in stream.");
        
        /* set number of neighbors value */
        outa = mxCreateDoubleMatrix(1, 1, mxREAL);
        if (outa == NULL)
            mexErrMsgTxt("Error creating neighbor list.");
        outd = (double *) mxGetPr(outa);
        if (outd == NULL)
            mexErrMsgTxt("Error allocating memory for neighbor list.");
        *outd = (double) nn;
        mxSetCell(*plhs, vn, outa);
        
        /* create double matrix */
        outa = mxCreateDoubleMatrix(1, nn, mxREAL);
        if (outa == NULL)
            mexErrMsgTxt("Error creating neighbor list.");
        outd = (double *) mxGetPr(outa);
        if (outd == NULL)
            mexErrMsgTxt("Error allocating memory for neighbor list.");
        
        /* copy numbers */
        for (; nn > 0; --nn)
            *outd++ = (double) (*ind++ + 1);
        
        /* set into cell */
        mxSetCell(*plhs, nv + vn++, outa);
    }

    /* set skip number */
    outd = (double *) mxGetPr(plhs[1]);
    *outd = (double) skip;
    
}
